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ABSTRACT 


Rapid flares from blazars in very high energy (VHE) y-rays challenge the common understanding of jets 
of active galactic nuclei (AGNs). The same population of ultra-relativistic electrons is often thought to be 
responsible for both X-ray and VHE emission. We thus systematically searched for X-ray flares at sub-hour 
timescales of TeV blazars in the entire Rossi X-ray Timing Explorer archival database. We found rapid flares 
from PKS 2005-489 and S5 0716+714, and a candidate rapid flare from 1ES 1101-232. In particular, the 
characteristic rise timescale of PKS 2005-489 is less than half a minute, which, to our knowledge, is the 
shortest among known AGN flares at any wavelengths. The timescales of these rapid flares indicate that the 
size of the central supermassive black hole is not a hard lower limit on the physical size of the emission 
region of the flare. PKS 2005-489 shows possible hard lags in its flare, which could be attributed to particle 
acceleration (injection); its flaring component has the hardest spectrum when it first appears. For all flares, the 
flaring components show similar hard spectra with I’ = 1.7— 1.9, and we estimate the magnetic field strength 
B ~ 0.1-1.0 G by assuming synchrotron cooling. These flares could be caused by inhomogeneity of the jets. 
Models that can only produce rapid y-ray flares but little synchrotron activity are less favorable. 


Subject headings: galaxies: active — galaxies: jets — X-rays: galaxies — gamma rays: galaxies 


1. INTRODUCTION 


Blazars, including BL Lac objects and flat-spectrum ra- 
dio quasars (FSRQs), are a special class of radio-loud ac- 
tive galactic nuclei (AGNs) that have one of their relativis- 
tic jets pointing very close to our line of sight (e.g., 
Padovani}1995). FSRQs have luminous broad emission lines 
that are weak or absent in BL Lac objects. Due to Doppler 
boosting, the emission of a blazar is usually dominated by 
the jet whose spectral energy distribution (SED) shows two 
broad humps that smoothly extend from radio to y-rays. The 
low-energy hump can extend from radio to soft X-rays. Ac- 
cording to the frequency of the first hump, BL Lac objects 
are further divided into low-frequency peaked BL Lac ob- 
jects (LBLs; Vpeak < 10'* Hz), intermediate-frequency peaked 
BL Lac objects (IBLs; 10!* < Vpeak < IO" Hz) and high- 
frequency peaked BL Lac objects (HBLs; Goen > 10!° Hz) 


E Abdo et al. 20703). The 
high-energy hump extends from hard X-rays to 7y-rays, even 
sometimes the very high energy (VHE) TeV band. Such VHE 
blazars, typically HBLs, are called TeV blazars. 

The low-energy hump is attributed to the synchrotron emis- 
sion of highly relativistic electrons gyrating in a magnetic 
field in the jet. The origin of the high-energy hump, however, 
is still debated. A popular explanation is inverse-Compton 
emission from the same population of relativistic electrons 
that produce the synchrotron emission. The seed photons of 
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the inverse-Compton scattering process could be local syn- 
chrotron photons (usually for BL Lacs) in the jet and/or ex- 
ternal photons from the central engine (usually for FSRQs) 
or the cosmic microwave background. These models are thus 
called leptonic models. The hadronic models, on the other 
hand, attribute y-ray emission to synchrotron emission of pro- 
tons (Mücke & Protheroe 2003} 
2015) or proton-induced cascades (Mannheim 
1998). 

The observed photon flux of blazars varies significantly 
across the electromagnetic spectrum on timescales from min- 


utes to years ee ee eee 
The origin of the variability 1s not well understood. Gener- 
ally, the variability is noise-like (e.g., 
Chatterjee et al.||2012), similar to the variability of radio- 
quiet AGNs (e.g., [Markowitz et al./2003). However, blazars 
are also known to have bursts that show flare-like structures 
(e. ns Sanaa 017 which may have recognizable pat- 
terns (Sasada et al.||2017). The outbursts can be explained 
by internal shocks of the jets when a new relativistic blob of 
plasma catches up with an old blob and accelerates particles 
to ultra-relativistic energies (e.g., EE Based 
on several similarities between the jet emission and corona- 
disk emission, the ultimate origin of the variability may still 
be accretion-rate fluctuations of the disk (McHardy|2008). 
The shortest variability timescale is a crucial parameter be- 
cause it serves as an independent constraint on the physical 
scale of the emission region (Tavecchio et al.|}1998), which 


cannot be easily provided by other observational measure- 
ments. Blazars are usually most variable at frequencies just 


above the two SED humps (e.g. [Uitich ot al [1777 [Madejski] 
& Sikora|2016), which usually fall in the hard X-ray and TeV 
bands in the case of TeV blazars (e.g., 
REN Tn particular, an 


increasing number of TeV blazars show y-ray flaring activ- 
ity on timescales from several to a few tens of minutes that 
are detected by ground-based Cherenkov telescopes, includ- 


ing both BL Lac objects (e.g., 1996} 


2007}|Albert et al.)/2007| 2013) and FSRQs 


(Aleksić 1 al.}2011 ve minute-scale variability in the TeV 
band (e.g. 2007) has strong EE for 
our EE of RGN jets (Begelman et al.!2008). X-ray 
and TeV emission may be directly related to E same high- 
energy tail of the relativistic electron population. Indeed, the 
lightcurves of HBL-type TeV blazars in the X-ray and TeV 
bands are usually correlated (e.g.,{Aleksié et al./2015b EK 

a extremely rapid X- 


niss et al./2015b). Attempts to searc 
ray variabi ‘Wek ave been made (e. z (ent es [Xue & Cuil 
Paliya et al. Pare BOTS} 


ute | SE Pryal et al 12018 e same source 
can have minute-scale variability ın both the X-ray band and 


TeV band (e. g., Mrk501, Xue & Cui}2005}|/Albert et al.|2007) 


However, “orphan” TeV e that ae no X-ray Gout 


parts are occasionally reported (e.g.,]/Krawczynski et al.|2004} 
Blazejowski et al.|2005}|Acciari et al.[2000}[Fraija et al |2015). 


Rapid TeV variability has germinated various models to ex- 
plain the small timescales. Most models involve some very 
compact regions moving in the rest frame of the jet. These 
compact regions could be “jets in a jet” that are either pro- 
duced by magnetic reconnection processes in a Poynting flux- 
dominated jet (Giannios et al./2009) or relativistic turbulence 
in the jet (Narayan & Piran|2012). The minijets-in-a-jet model 
can consistently produce the statistical properties of blazar 
flux (Biteau & Giebels|2012). There are also models involv- 
ing a red giant star being stripped of its envelope by the jet 
(Barkov et al.|2012) and models involving beams of magneto- 
centrifugally accelerated electrons occasionally pointing to- 
ward us (Ghiselin et al.[2009) 

Several well-studied TeV blazars show rich spectral behav- 
ior in X-rays, which may represent the general behavior of 
the synchrotron peak of all AGN jets. The X-ray spectra are 
usually curved (Massaro et al.!2004) and can only locally be 
fitted by a power-law. 7 e spectral at with flux can 


be complex (Zhang et al. AER (Cui|/2004). Generally, the 
spectrum hardens when the flux increases (eg. 
21006} Kue et al.)2006) 7 12009), but photon in 


dexes can saturate at higher fluxes (Xue & Cui]2005} [Giebels] 
et al.2007). The synchrotron peak usually moves to DE 


frequencies with increasing flux during outbursts (e.g.|P 
1998), but no correlation between the break energy and the 
Hux exists when a broken power law is adopted to fit the e 


ray neg OIT (Xue & Cuif2005} |Giebels et al./2007| 


2010). A cooling break in the Dee of emitting = 


a es SS ei explain these features (Wierzcholska & Wag- 
mer 20161, and some special particle acceleration processes 


may be involved (Madejski & Sikoral2016). There are also 
energy-dependent lags between the variations of different en- 
ergy bands. In some flares, soft bands lag behind hard bands 


(e.g., Zhang et al.|2002), while lags "7 me Opposite direction 
can also happen (e.g. [Ravasio et al 2008) Sato et al 2008). 
Hysteresis in the HR (hardness ak; —flux diagram 1s often 
used as a diagnostic of lags. Clockwise loops (e.g.,|Acciari 
e al.{2009) Kapanadze et al-2016} in the HR-fux plane are 
a sign of soft lags while counterclockwise loops (e.g. 
2009) are a sign of hard lags. The same source 
can exhibit both clockwise and counterclockwise loops; the 
observed patterns are further complicated by the superposi- 
tion of flares at different timescales (Cuil2004). The above 
knowledge of TeV blazars in the X-ray regime comes from 
studies focusing on timescales of hours to weeks. We will ex- 
tend this kind of analysis to much smaller timescales in this 
paper. 

The main goal of this paper is to search for X-ray flares 


TABLE 1 
TEV BLAZARS WITH > 50 RXTE/PCA OBSERVATIONS 

Name Z Type Numberof Exposure 
Observations Time (ks) 

3C 279 0.5362 FSRQ 1988 3198 

BL Lacertae 0.059 IBL 1387 2522 

Mrk 421 0.031 HBL 1190 2515 

PKS 1510-089 0.361 FSRQ 1334 2254 

PKS 2155-304 0.116 HBL 501 1107 

Mrk 501 0.034 HBL 499 886 

S5 0716+714 0.31 IBL 233 733 

H 1426+428 0.129 HBL 164 527 

PKS 2005-489 0.071 HBL 158 483 

3C 66A 0.41 IBL 99 373 

PKS 1424+240 vee HBL 64 347 

1ES 0229+200 0.14 HBL 205 295 

1ES 1959+650 0.048 HBL 147 272 

1ES 1101-232 0.186 HBL 99 211 

1ES 2344+514 0.044 HBL 53 134 


NOTE. — We list above all the TeV blazars with > 50 PCA observations 
that add up to > 130 ks exposure time. The remaining unlisted TeV blazars 
have < 50 PCA observations, which are PG 1553+113, 1ES 1218+304, 
MAGIC J2001+435, 1ES 0806+524, 1ES 0647+250, RGB J0152+017, 
1ES 0414+009, W Comae, LES 1727+502, Mrk 180, PKS 0447-439, 
RGB J0710+591, PKS 0548-322, AP Librae, LES 1741+196, and 

H 2356-309. See for the full list of known TeV 
blazars and their redshifts and classifications. 


at sub-hour timescales from TeV blazars in the entire Rossi 
X-ray Timing Explorer (RXTE) archival database. We use 
data from the narrow-field pointing instrument Proportional 
Counter Array (PCA) onboard RXTE, covering a nominal en- 
ergy range of 2—60 keV. The RXTE satellite was launched 
in December 1995 and ceased science operation in January 
2012. During its lifetime, it accumulated more than ~ 16 Ms 
of exposure time on TeV blazars in hard X-rays, surpassing 
any other X-ray observatory. We describe data reduction and 
the searching results in Section|2| Most of our following anal- 
ysis is based on an assumption that the observed photons are 
from a flaring component and an underlying constant/slowly- 
varying component, possibly from two separated sites. We 
describe lightcurve-model fitting and spectral-model fitting in 
Section [3] and Section [4] respectivi. We discuss the impli- 
cations of our findings in Section [5] and summarize them in 
Section |6| In the following, we use the ACDM model, with 


Ho = 67.7 km/s/Mpc and Qm = 0.307 (Planck Collaboration 


2016). 


2. DATA REDUCTION AND SEARCHING FOR FAST FLARES 


We retrieved all the archival RXTE/PCA observationd’] of 
TeV blazars (see Table DI The total number of observations 
is ~ 8400, and the total exposure time is ~ 16 Ms. We used 
the Standard 2 mode data, which have a time resolution of 
16 s._The data were reduced following the standard proce- 
dure|* We created filter files and good time intervals (GTIs) 
for each observation according to the suggested screening cri- 
terid?| for faint sources. Background data were then simulated 
using the appropriate model "| We applied the GTIs to both 


7 The data were downloaded from ‘http://heasarc.gsfc.nasa.gov/cgi- 


bin/W SE d 


background model file 
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FIG. 1.— Rapid X-ray flares of TeV blazars in bins of 16s. The red dashed 
curves are the weighted least square models (Eq. HI The vertical dashed 
lines are used to separate the flaring phase from the quiescent phase in the 


light curve if possible. The event of Mrk 501 was first reported We 
(2005), which we reanalyzed more quantitatively in this paper. e count 
rates for the y-axes are normalized to one PCU. The corresponding ObsIDs 
of the three events from top to bottom are 20342-03-01-01, 95377-01-91-00, 


and 30249-01-01-02. On top of each panel, the date (Modified Julian Date) 
when the flaring observation started (Oe, set as t = 0 s) is annotated. 


observational data and simulated background and extracted 
lightcurves in initial 16s bins from channels that correspond to 
~2-20 keV. The lightcurves of net count rates were calculated 
using 1cmath in the HEASoft (v6.19) package. Since only 
PCU2 among the five proportional counter units (PCUs) of 
PCA is almost always in operation, we extracted lightcurves 
from PCU2 for flaring event selection. We visually inspected 
every lightcurve to select events in individual observations 
that contain a complete or nearly complete sub-hour flaring 
profile. Specifically, we require the flare to have apparent rise 
and fall; we also require the existence of a plateau either be- 
fore the rise and/or after the fall to assess the completeness of 
the flare and the level of the background component)’ | Note 
that RXTE is in low Earth orbit, and thus an uninterrupted 
lightcurve is usually less than ~ 50 minutes due to the Earth’s 
occultation or passage through the South Atlantic Anomaly, 
etc., which limits the timescales of the events investigated. 
After identifying fast flaring events, we extracted lightcurves 
and spectra from all the PCUs available during that observa- 
tion to achieve higher signal-to-noise ratio (S/N) for further 
analysis. 


'l We do not adopt quantitative criteria to select the events, automatically, 
because quantitative criteria are inescapable of subjective tweak and visual 
inspection is almost always necessary. 
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From the complete RXTE/PCA database, two new fast flar- 
ing events were found. They belong to PKS 2005-489 and 
S5 0716+714. The lightcurves of these two events are shown 
in Fig.[I] We also plot an event of Mrk 501 (reported in|Xue &| 
2005). A fast-flaring candidate event of 1ES 1101-232 is 
reported in Appendix [A] which has relatively low credibility 
because of limited S/N. These flaring observations generally 
lack simultaneous observations in other wavebands. 

We checked for potential contamination by soft electron 
flares that were not screened out by the criteria “ELEC- 
TRON2.LE.0.1” in data cleaning. The contemporaneous 
Electron2s of each event were well below 0.1 and did not 
show any apparent electron flaring activity that may be re- 
sponsible for the X-ray flares. The longitudes and latitudes of 
the satellite at the onsets of the X-ray flares did not cluster in 
the anomalous high background region (cf. Fig. 7 of 


2005). Further support for the genuineness of the flaring 


events comes from their lightcurve and spectral features ex- 
plored below. They behave like well-known X-ray flares of 
TeV blazars, only at much smaller timescales. In conclusion, 
we did not find any sign of contamination of soft electrons or 
any other known sources for all the flaring events. However, 
the possibility of an unrelated X-ray transient in the field of 
view still cannot be ruled out entirely since the PCA lacks the 
capacity of imaging. 


3. LIGHTCURVE FITTING 


Lightcurves in additional energy bands (see the second col- 
umn of Table [2] and Figs. are_extracted according to 
the energy-channel conversion "d Thanks to the high 
data quality, variability at timescales down to the time res- 
olution (16s) of the lightcurves (see Fig. is seen. Be- 
cause of the variety of data quality and gain epochs (Jahoda 


2006), we do not have uniform definitions for different 


bands for all three sources. We have four different bands for 
PKS 2005-489, while we have three bands for Mrk 501 and 
S5 0716+ 714; the full band of PKS 2005-489 and Mrk 501 
is 1.94—20.30 keV, while the full band of S5 0716+714 is 
2.06—10.11 keV. 


3.1. Method of Fitting 
We fitted the lightcurves with a constant flux plus an ex- 


ponentially rising and decaying flare following (Abdo et al. 
(2010b}!> 


=f t—t0 


i -1 
FO=R+R (eF +e) , (1) 


where F, represents the constant flux level underlying the 
flare, and 7, and 74 are the characteristic rising and decay- 
D D 14 D D . o 

ing timescaleq!4| of the flare. tọ indicates the transition from 
rising to falling, and the count rate actually peaks at 


EEN (2) 


Tr + Td Tr 


which equals fo only when the flare is symmetrical (7, = Ta). 
Therefore we define the amplitude of the flare as F,/F,, i.e., 
the count rate at fo (Fp) over the constant level, instead of 


12 See the table at[hup//heasare-gsfe nasa;gov/docs/xte/e-c_table-html 
13 Tn addition to the frequently-used Eq. T some similar analytical expres- 

sions have been used to describe the flare profiles of blazars; see[Albert et al 
+ The doubling and ha x 


ving timescales are 7; X In2 and Tq x In2. 


TABLE 2 
FITTING RESULTS FOR THE FAST FLARING EVENTS 
Band Energy Fe Fo to ge Ta tp Fy/Fe E x2, /dof 
(keV) (c/s) (c/s) (s) (s) (s) (s) 
PKS 2005-489 
Soft 1.94-5.47 2.0819% 3.814058 30318 234 IER 2437 2.3 0.75  0.957/166 
Medium 5.47-10.11 0.9819%  3.45+059 31810 29+% 141414 3567 3.2 0.66  1.05/166 
Hard 10.11-20.30 0.337705 e, 324133 Boke ‘105% 353777 53 0.56  0.895/166 
Full 1.94-20.30 3.41904 94904 31118 254 143+! 2487 2.8 0.70  1.18/166 
S5 0716+714 
Soft 2.06-5.31 0.197793 2.044938 (äi 222" 82546 Lamm 75 0.58  1.48/23® 
Hard 5.31-10.11  0.00t9:95 1.672931 1216+% 161s Issey. 1551419 >17® 0.82 16203” 
Full  2.06—10.11 OI 3.95+058 1258+37 17840 1158+320 1547+% 17 0.73 1.15/159®) 


Mrk 501 first substructure 


Soft 1.94-5.82 5.60720) 3.321902 1918133 9977 11072 1924+14 Län 0.05 0.99/151 
Hard  5.82-20.30 5.35400) 3.297028 200148, 1471s 33165 196178, 1.38 0.64 1.09/151 
Full 1.94-20.30  10.96+9:06 6.970% 1961418 (ol 67i 19361? 1.33 -0.28 1.23/151 


Mrk 501 second substructure 


Soft - - 17803 220659 0 — 160%] - - 0.60 - 
Hard - - 73, Se, ww ww - - 0.82 - 
Full : : Be ie oe SES S S 0.67 S 


"Note our fitting cannot constrain the hard-band quiescent flux level of S5 0716+714, so the hard-band amplitude is unmeasurable. We thus set the full-band 


amplitude as the lower limit of that of the hard band. 


>We fitted the model to the soft- and hard-band lightcurves in bins of 96s, while the full-band lightcurve is in bins of 16s. 
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FIG. 2.— Top three panels: The rapid X-ray flare of PKS 2005-489 in three 
energy bands in bins of 16s. The red solid curves are weighted least square 
models; the vertical blue dashed lines indicate tp, whose 1o uncertainties are 
shown as the shaded blue regions. Bottom two panels: Hardness ratios (hard 
band to soft band and medium band to soft band); the red solid curves are 
hardness ratios calculated from the weighted least square models above. The 
data points are in bins of 32s. 


Fo/F-. The symmetry of a flare is described by 


Td — Tr 


Tat Te 


(3) 


whose value is in the range of [—1, 1]. € =—1 (= 1) represents 
completely right (left) asymmetric profiles with a zero falling 
(rising) timescale; € = O indicates a symmetric flare. 

The uncertainties of t, and F, have to be propagated from 
the errors of other parameters. We adopted an MCMC 
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FIG. 3.— Soft- and hard-band lightcurve fitting of S5 0716+714. The 
lightcurves are in bins of 96s. The vertical dashed lines indicate tp, whose lo 
uncertainties are shown as the shaded blue regions. 
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FIG. 4.— Soft- and hard-band lightcurve fitting of Mrk 501. The lightcurves 
are in bins of 16s. The vertical dashed lines indicate the peak of the main flare 
tp, whose lo uncertainties are shown as the blue shaded regions. 


(Markov Chain Monte Carlo) algorithm to fit the lightcurves, 
which returns reliable probability intervals of timescales and 


amplitudes by sampling from their posterior distributions. We 
first performed weighted least squares fitting using a numer- 
ical minimizer to obtain the best estimates of F}, Fo, to, Tr, 
and Ta. Starting from these initial values, we took 1000 ran- 
dom walk steps in parameter space. The samples of t, and Fp 
were calculated according to Eq. and Eq.}1} Note that we 
added a second flaring component to the model in the fitting 
of Mrk 501 (see the bottom panel of Fig. |). 


3.2. Correcting Lightcurve Error Bars 


The initial fits have reduced Chi-square values (x2 = 
vi /dof) in the range of 0.48—-0.80, which indicates that the 
assigned error bars are larger than true statistical fluctuations. 
Indeed, the standard RXTE/PCA data reduction pipeline over- 
estimates the lightcurve errors (Nandra et al 20000. The error 
estimation of the net lightcurves is propagated from the er- 
ror estimation of the observed lightcurves and the simulated 
background lightcurves, of which the latter is too smooth to 
be described by the assumed Poisson statistics. We decided to 


correct the error estimation using oj, = Oe TF E Gite, where 


0 <k* <1. The correction factor k? can be determined by 
forcing the excess varianc d" of the quiescent parts in the top 
and bottom panels of Fig.|1|to be zero. The resulting cor- 


rection factors of different segments at different energy bands 
span from —0.13 to 0.71. We decided to fix k? =0 as in Nand 
er (2000). We have ignored errors on the background in the 
ightcurve analysis below, unless otherwise stated. The fitting 
results after correcting the error bars are tabulated in Table 
Note that we still report the least square results in Table [2)as 
the estimation of each parameter, but the Io intervals are de- 


rived from MCMC fitting. The reduced Chi-square values of 
most fits are around 1. 


3.3. Lightcurve Fitting Results 


Every lightcurve shows some flare-like structure above a 
constant “background” flux level, which actually varies on 
longer timescales. From the flux levels of the constant compo- 
nents, the events occur when the sources are in relatively high 
states, but they do not always coincide with the periods with 
the highest flux levels. For example, the RXTE/PCA count 
rates of PKS 2005-489 can be 10 times the constant flux level 
here as found about one and a half years later 
1999). 

The variation amplitude is higher in harder bands, which 
suggests that the flaring component has a harder spectrum 
than the corresponding constant component. The mixed spec- 
tra become harder when the flux rises and the flaring compo- 
nent becomes more prominent. We discuss the spectral vari- 
ability of PKS 2005-489 in detail in Section We also fit 
the spectra of PKS 2005-489, S5 0716+714, and Mrk 501 
and confirm that the flaring components have harder spectra 
in Section [4] 

There is a trend of rising timescales being shorter and de- 
caying timescales being longer at softer energies. This trend 
is obvious in PKS 2005-489 and the first flare substructure 
of Mrk 501 (see Table [2} see also Mrk 421 in Appendix B). 
This suggests that the variability is caused by electron accel- 
eration and cooling (eg. |Fraija et al.|/2017). As a conse- 
quence, the flares are more right asymmetric in harder bands 
(i.e., smaller € values). S5 0716+714 and the second flare 


15 Excess variance is the variance after subtracting the mean square error 


(e.g., 1997\| Vaughan et al./2003). 
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substructure of Mrk 501 do not follow the patterns, although 
the error bars prevent any solid conclusions. The differ- 
ence of timescales between energy bands is not obvious for 
S5 0716+714, mainly because the decay of the flare was not 
completely sampled and 7g cannot be constrained well. More- 
over, only one PCU was operating during the observation, so 
the S/N is low. 

The vertical lines in Fig. [2|and Fig. |4] suggest hard lags in 
the variability of PKS 2005-489 and Mrk 501, which means 
that the variation of hard photons lags that of soft photons. 
The suggested lag of PKS 2005-489 is not confirmed by 
cross-correlation function (CCF), presumably due to the lim- 
ited time-resolution (see discussion in Section [3.4). Again, 
the soft-band and hard-band peaks of S5 0716+ 714 do not 
show an obvious difference, putatively due to the larger error 


bars (see Fig. D. 


3.4. Spectral variability of PKS 2005-489 


The flare of PKS 2005-489 has the shortest rising 
timescales and its data are of the highest quality, so we inves- 
tigated its spectral variability in detail. We plot the hardness 
ratios! (HRs) in the bottom two panels of Fe UU which show 
hardening that corresponds to the flare. We also plot a HR- 
flux diagram of this event in panel (a) of Fig.[5] which shows a 
“harder when brighter” trend and hysteresis. The loop begins 
with clockwise motion, and then follows a counterclockwise 
direction. Below, we argue that the clockwise trend at the 
beginning is due to the superposition of two spectral compo- 
nents. 

The apparent two-component nature of the lightcurves sug- 
gests that the spectral variability is partially caused by a 
change of the relative fraction of the two components. In 
other words, if the flaring component has different hardness 
ratio from that of the constant component, then even if nei- 
ther of the spectra changes over time, the observed overall 
hardness ratio will still change due to the flux variation of the 
flaring component (see 2014) 2015). 
However, the hysteresis loop in panel (a) suggests the spec- 
trum of the flaring component is intrinsically variable; other- 
wise the track in the HR us plane while the flux is rising 
will be identical to the track while the flux is declining, in- 
stead of forming a loop. To discriminate the effects caused 
by the mixing of different components and the flare’s intrinsic 
spectral variability, we further subtract the constant compo- 
nent F, from the lightcurve of each band, where F, is from 
Table We calculated the hardness ratios from the result- 
ing flare-only lightcurves in panels (c)(d)(e) and HR—flux di- 
agram in panel (b). The plots suggest that the flare emerges 
with a hard spectrum, softens gradually as the flux rises, then 
hardens near the highest flux, and finally softens as the flare 
fades away. Only a counterclockwise loop is apparent in the 
HR-flux plane for the flare-only lightcurve. Based on the ev- 
idence above, we know that the clockwise loop in panel (a) 
could be caused by the sudden emergence of a hard flaring 
component and the spectral variability soon follows the spec- 
tral variability of this flaring component due to its increas- 
ing dominance. The flaring component itself has complex 
spectral variability and the counterclockwise loop in panel (b) 
suggests the existence of a hard time lag, albeit being sub- 
ject to large uncertainties. Note that the overall oblique “8” 


16 We define hardness ratio as the count rate of hard band over that of soft 
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FIG. 5.— Panels (a)(b): HR-flux diagram of PKS 2005-489 before and after subtracting the constant component F}. Blue arrows indicate the time sequence. 
The effective photon index I’ shown as the y-axis on the right hand side is obtained using the response files of the PCA, a Galactic absorbed power-law model, 
and a range of assumed photon indices. Panels (c)(d)(e): Three hardness ratios of the flare versus time after subtracting the constant component Fc; the vertical 
dashed lines indicate the time of full band peak (tp). We only considered the ~ 400s (from ~ 200s to ~ 600s) segment which contains nearly the whole flare. 
The data points are in bins of 32s. The dashed curves are calculated from analytical models (Goal 


shape before subtracting the constant component in panel (a) 
is reminiscent of a flare of Mrk 421 which lasted ~ 60 ks 
and was reported EE (2010) using Suzaku data. 
We reanalyzed the Suzaku data of Mrk 421 and performed 
the same spectral variability analysis as for PKS 2005-489 
above in ee a We reproduced the same results as for 
PKS 2005-489 in the flare of Mrk 421 and thereby strength- 
ened our conclusions. The similar behavior of spectral vari- 
ability, seen both in the extremely rapid flare of PKS 2005- 
489 and in the long-duration flare of Mrk 421, is reminiscent 
of the scale invariant nature of X-ray flares from TeV blazars 


(Cuil2004}/Xue & Cuil2005). 


We also calculated the CCF (Edelson & Krolik| 1988} 
1999) of the soft band and hard band flare-only lightcurves, 


but did not find an obvious time lag. Any time lag in this fast 
flaring event of PKS 2005-489 could be intrinsically small 
compared with the time resolution of the observation (16s). 
The bandpass of our data is narrow, spanning about one or- 
der of magnitude, so any energy-dependent lags may not be 


significant. The size of lags may also be positively correlated 
to the duration of the flares EE olf which means 
that flares of shorter duration have smaller lags. Note that in 
Table|2|the differences of fo, tp, and 7, between different bands 
are small compared with relatively large error bars. 


4. SPECTRAL FITTING 


Motivated by the two-component model in lightcurve fit- 
ting, we carried out spectral fitting in a similar way. By divid- 
ing the observation into flaring and quiescent phases accord- 
ing to the lightcurves, we separately extracted and jointly fit- 
ted the spectra of the two phases. The spectra were extracted 
from the top layers of the operating PCUs. The spectra and 
models are shown in bg Gomm the E Fg representation and the 
photon indexes are tabulated in Table |3| We used channels 
that correspond to 3—20 keV for PKS 2005—489 and Mrk 501 
and 3—10 keV for S5 0716+714, respectively. In addition 
to Galactic absorption, we used a power law and sum of two 
power laws in XSPEC (Vv 12.9.0) to fit the spectra of the qui- 
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FIG. 6.— In the spectral fittings (the top row and bottom-left panel), the Galactic absorption is fixed to Ny = 5.08 x 10% cm, 3.81 x 107° cm~? for 
PKS 2005-489, S5 0716+714 and 1.56 x 102° cm? for Mrk 501 (Kalberla et al.!2005). Note that the average spectra of the 
flaring components are shown in green dash-dotted lines. In the bottom-right panel, we show the observations that were taken before the flare (pre) and after the 
flare (post) of S5 0716+714, whose ObsIDs are 95377-01-90-00 and 95377-01-92-00. Note that the vertical dashed lines represent gaps of 3—4 days between 
adjacent observations. We fit the pre-flare and post-flare lightcurves with constant fluxes and show the results by horizontal lines. Shaded regions indicate lo 


uncertainties. 


TABLE 3 
SPECTRAL FITTING RESULTS OF THE FAST FLARING EVENTS 
Source IR I% F ies Fiare ® Fiare onl (a) x, /dof Lx® Ed 
PKS 2005-489 2.73+0.05 1.7140.18 29.5 42.4 12.9 0.88/38 3.86x IO 2.08 x 10% 
Mrk 501 2.01 +0.02 1.86 +0.20 115.7 130.0 14.3 0.76/38 3.23 x 10% 5.13 x 10% 
S5 0716+714 2.48 +0.27 1.89+0.15 5.7 (pre) 12.3 12.2 0.32/16 1.90x10* 1.06 x IO 
3.8 (post) 1.25 x 10% 


Fluxes are in units of 1071? ergs cm~? s7! ; the energy range is 3—20 keV for PKS 2005-489 and Mrk 501, but 3-10 keV for S5 0716+714. 


>X-ray luminosity of the constant component luminosity in units of ergs s71. 
“Energy of the flaring component in units of ergs. 


escent and flaring phases, respectively. We tied the power law 
in the quiescent phase to one of two power laws in the flaring 
phase, for both normalization and photon index (see the leg- 
ends of Ps loi, We also calculated the X-ray luminosity of the 
constant component and the total energy of the flaring compo- 
nent in the X-ray band (namely average flare luminosity times 
the length of the defined flare phase) in eet Note that the 
spectra should suffer from overestimation of error bars due 
to the same reason discussed in Section [3.2] but we are not 
able to correct the error estimation of the spectra as we did 
in the lightcurve fitting. One consequence of the overestima- 


tion of error bars is that the parameters in the spectral fitting 
have large confidence intervals. Another consequence is that 
we might be able to fit the spectra with many models. For 
example, we can fit the spectra of the flaring phase, quies- 
cent phase, and whole observation with simple power laws, 
and the Chi-square statistics can still be acceptable. We think 
that describing the spectra as a combination of different power 
laws is more physically appropriate and is consistent with the 
lightcurve fitting. 

Note that due to the limited length of the observation, 
we have neither a complete flaring nor a quiescent phase of 
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S5 0716+714. Therefore we extracted the spectra of the two 
observations just before and after the flare, which are referred 
to as pre-flare and post-flare, respectively. The gaps between 
the observations are 3—4 days. The three observations are 
plotted in the bottom-right panel of Fig. o where we fitted 
each of the pre-flare and post-flare lightcurves with a constant 
flux. The spectra are jointly fitted with a tied “constant” com- 
ponent photon index, but the normalizations of the “constant” 
component are free parameters in this case. In other words, 
the constant component varies on timescales of a few days. 
If we leave the normalization of the constant component un- 
derlying the flare (factor “A” in the legend of top-right panel) 
free, the value of “A” in the fitting results is negligible but its 
error bars are large. Additionally, we are not able to constrain 
the photon index of the flaring component (l'2) very well in 
this case. Thus we fixed A = 0.6D, according to the lightcurve 
fitting results in the bottom-right panel of Fig. [6] The choice 
of 0.6 or another reasonable value does not affect the result of 
Ts <T). 

The above fitting shows that for all three flares, the flaring 
components have harder photon indexes than the correspond- 
ing constant components (see Table (3). The difference is not 
as apparent for Mrk 501 because of the relatively small varia- 
tion amplitude (e, relatively weak flaring component). The 
photon index of the constant component ranges from 2.0 to 
2.7, indicating that we are observing different declining parts 
of the synchrotron hump for different sources. The photon 
indexes of the flaring components, on the other hand, lie in a 
narrower range around 1.8. Every estimated average flux (av- 
eraged over the length of the flare phase we defined) of the 
events is ~ 107!! ergs cm? el (see Table|3). Using the “effi- 
ciency limit” of compact sources (e.g.,/Fabian| 1979} 
al.|1999), AL/At < 2 x 10*!no.1 erg s™!, the matter to energy 
conversion efficiency 7 is greater than 1 for PKS 2005-489 
and S5 0716+714, suggesting the existence of significant 
boosting of the emission (e.g., [Remillard et al./1991). Due 
to the small flare amplitude, the efficiency of Mrk 501 does 
not exceed 1. 


5. DISCUSSION 
5.1. The Flaring and Constant Components 


We interpret all the lightcurves that we analyzed in detail 
(including Mrk 421 in Appendix |B) as a superposition of a 
constant component and a flaring component (e.g.,|Fraija et al. 
(2017). The two-component model depicts a simpler scenario 
than a single component undergoing an outburst. For the case 
of a single emitting region, the variation in the HR—flux plane 
is more complex. Furthermore, there must have been a sudden 
enhancement of some key physical quantity and this quantity 
has to fall back later to its value preceding the burst. 


5.2. X-ray Radiation Process of A3 0716+714 


S5 0716+714 is an IBL (Ackermann et a 201 In with a syn- 
chrotron hump peaking at optical wavelengths (Anderhub et 
REES It has frequent intra-day variability and hysteresis 

oops 1n the color-magnitude plane in the optical band (eg, 
Man et al.{2016). (2015) reported two rapid flares 

95 0716+714 in the X-ray band with low significance. The 
X-ray emission covered by the RXTE bandpass could be a 
mixture of both synchrotron and inverse-Compton emission 
(Wierzcholska & Siejkowski[2016). As such, S5 0716+714 
will not necessarily behave the same way as HBLs in X- 
rays. Depending on the state (low or high) of the source 
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FIG. 7.— Variability timescale versus supermassive black hole mass. We 
use the rising timescales in the galaxy’s frame, 7; /(1+z), as the values on the 
y-axis. For convenience, we annotate on the right y-axis the corresponding 
sizes, T-c/(1+z), in units of cm and gravitational radius for Mgy = IO AM. 
We also calculate light-crossing times for Kerr black holes, tj. = 2GMpy / C= 
2x 10° (Men / 10°Mo)s. The lines represent light-crossing times after con- 
sidering the relativistic Doppler effect, 4. /6. References of vee Casco 
(2008), (2) (2002), (3) (2009), (4) 
when the observation was made, the dominant X-ray radi- 
ation process may change accordingly. The spectrum may 
steepen when S5 0716+714 brightens since the synchrotron 
emission in the soft X-rays becomes increasingly important 
(Giommi et al.J1999). However, the synchrotron component 
itself hardens due to the shifting of the synchrotron peak 
EE et al. [2006] [Zhang] 2010). The steep photon index and 
the high flux levels (3—10 keV) of pre-flare and post-flare ob- 
servations are consistent with synchrotron emission. The flar- 
ing component has a slightly inverted spectrum, but the even 
higher flux level and short timescales still support an origin as 
synchrotron emission. Therefore the flare we observed from 
S5 0716+714 does not seem to have a different radiation pro- 
cess from HBLs. 


5.3. Constraints on the Emission Region 


The upper limit on the physical scale of the flaring region 


can be given by: 
A tA 
ml ( ect 4 
(5) rg ee H 


where ô is the Doppler boosting factor, and tfare is the ob- 
served variability timescale. Since the flares are asymmet- 
ric, the size of the emitting region is reflected by the rising 
timescale 7, (Zhang et al.|/2002). Searching the literature, 
we find black-hole mass estimates for three out of the four 
sources (except LES 1101-232). We plot the rising timescale 
in the galaxy’s frame against black-hole mass in Fig.|7| Note 
that the typical Doppler factor of BL Lac objects is 10—20 
(20-30 for FSRQs) (e.g. EE The size of 
the black hole is often thought to be a natural lower limit on 
the physical scale of the emission region. We also plot this 
lower limit as a function of black-hole mass assuming dif- 
ferent Doppler factors. The variability timescales should lie 
above the lines, if the size of the black hole is a hard lower 
limit. The rising timescale of PKS 2005-489 is, as far as we 


Se CtflareO as 
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know, shorter than any variability timescale of AGNs at any 


wavelength ever reported (e.g., Remillard et a 11901 
1997} |Gallo et al.|2004}|/Xue & Cu1|2005}/Aharonian et 


al.||2007; (Albert et al.|/2007}| |Aleksi¢ et al.|/2014) 
. Fig.|7|/shows that the Doppler factor of PKS 2005-489 


has to be larger than several hundred, which appears unrealis- 
tic, in order to support the idea that the black hole sets a lower 
limit on the physical size of the flaring region. So far, there 
is no correlation found between the observed minimum vari- 
ability a and the black-hole mass 
E Babié/2015). The events of PKS 2005-489 and Mr 
(Appendix |B) dréien remarkably similar spectral Suton 
supporting the same process driving the spectral variability. 
However the timescales of the two events are different by two 
orders of magnitude, reflecting an intrinsic difference between 
the two jets. Noticeably, the famous Mrk 421 has many more 
X-ray observations than PKS 2005-489; however, the vari- 
ability timescale of Mrk 421 has never been found to be as 
short as the PKS 2005-489 event 
Paliya et al.[2015). A lower limit on the variability timescale 
probably does exist, and probably is not set by the central su- 
permassive black hole. 

The synchrotron cooling time of emitting electrons is given 
bY Lol Si 6rmec/ oryB? (Rybicki & Lightman (Rybicki & Lightman & Lightman| 1979) |1979} ). The 


Ges photon energy at the synchrotron peak is given b 
= hv = Beh/4rm,c)dy7B (Rybicki & Lightman & Rybicki & Lightman|1979). 
Canna the above two equations, we have 


tooo = 3.04 x IO BA? 6 AEN? s (5) 


where Ep is in units of observed keV (Zhang et al.|/2002). 


If we take PKS 2005-489 as an example, and adopt Lat = 
Ta ~ 143 s, Ep = 10 keV and ô = 30, we can have B ~ 1.2 G. 
The estiniited magnetic field of Mrk 501 is similar to that of 
PKS 2005-489, but the magnetic field of S5 0716+714 is 4 
times weaker. Radio-loud AGNSs are potential accelerators of 
cosmic rays. We estimated the maximum energy of protons, 
if they can be accelerated in the same region as the electrons, 
as Emax = eBR ~(2-9)x10!° eV. 


5.4. The Rarity of Extremely Rapid Flares of Blazars 


It is surprising to find such rapid variability of 
PKS 2005-489, because it was not found to be variable on 


timescales less than a day 1999} [Rector & Perl- 


[man][2003} HESS Collaboration et al.]/2010] [2011). Out 
of the ~160 PCA pointings_on PKS 2005-489, only two 


have ES eXCess varianed!”| U| including the one shown in 
Fig.|1{a). 

SC although the list of blazars that show extremely 
rapid flares is growing, these events are rare (e.g., 
1986) and usually unexpected. The rate of occurrence 
of sub-hour flares is about once per 4 Ms in the RXTE/PCA 
database. Most of the sources have only one such event re- 
ported, either in y-ray or X-ray (e.g.,|Gaidos et al.|1996}|Aha-| 
2007} and this paper). These rapid flares do not 
seem to be the extremely short cases from a continuous dis- 
tribution of the timescale of flares (e.g.,{Li et al./2017}/S 
et al.[2017). It remains unknown whether such events exist in 
all wavebands, so we do not know the total energy output of 
the flares. The biggest challenge is to coordinate multiple in- 
struments to target the same source simultaneously and hope 


17 We do not correct the error bars of the light curves here, but the light 
curves are consistent with being flat under visual inspection. 


rare unpredictable flaring events happen. 


5.5. Particle Acceleration 


Assuming the fast-cooling regime, the photon index of the 
flaring component OU =1.7—1.9) indicates that the accelerated 
electrons have an effective energy spectral index (p = 2I'—2) 
in the range of 1.4—1.8, which can be achieved by relativis- 
tic magnetic reconnection (e.g., 2014). Indeed, the 
SED of PKS 2005-489 shows very low Compton dominance 


2011), which indicates 
high magnetization. The results of SED modeling (e.g., 
derhub et al.]2009} 
et al./2011)|Aleksi¢ et al./2015a), which usually assume one 
single homogeneous emitting region, have low-strength mag- 
netic fields (B ~ 0.01—0.1 G) compared with the estimation 
of Section|5.3}and below equipartition. This may require that 
either the jet is structured (Ghisellini et al./2005) or only the 
regions that are responsible for the fast flares we observed 
have high magnetization. 

The direction of the loops in the HR us plane is thought 
to be determined by the competition between the accel- 
eration/ejection timescale, cooling timescale, and escape 
timescale (Kirk et al.|/1998). In practice, it is actually de- 
termined by which part of the synchrotron spectrum we are 


observing. The distinctive pattern in the spectral variation of 
the flares of PKS 2005-489 and Mrk 421 is not predicted by 


time-dependent homogeneous one-zone models (e.g., 
variabil- 


1998} |\Chiaberge & Ghisellini}1999). The spectral 
ity pattern in Figs. |5)and|B2| could possibly be produced by 


time-dependent inhomogeneous blazar models (e.g., 
2010). 


6. SUMMARY 


We searched the entire RXTE archival database for rapid 
X-ray flares of TeV blazars that last less than one hour. 
We investigated the temporal and spectral properties of the 
fast flares discovered under a two-component assumption. 
Our analysis also includes an X-ray flare of Mrk 421 us- 
ing Suzaku data that has similar spectral variability to that 
of PKS 2005-489 in Appendix |B] Our main findings are as 
follows: 

1. We discovered two new fast X-ray flares from 
PKS 2005-489 and S5 0716+714 and a candidate flare from 
IES 1101-232. The event of PKS 2005-489 shows, as far 
as we know, the most rapid variation of AGNs that has been 
observed at any wavelength. The extremely small timescale 
(7, < 30 s) defies the size that corresponds to the light- 
crossing time of the supermassive black hole as a lower limit 
on the size of the flaring region. 

2. The flares are usually superimposed on a 
constant/slowly-varying component. The flaring com- 
ponent generally has a harder X-ray spectrum (I = 1.7—1.9) 
than the constant component (1° = 2.0-2.7). 

3. The higher data quality of the X-ray observations can 
provide more detail than y-ray observations. The flare-only 
component shows a counterclockwise pattern in the HR—flux 
diagram, providing a sign of likely hard lags. This component 
also has the hardest spectrum right at its appearance. 


RXTE ceased science operation in January 2012, but X-ray 
observatories like Chandra, XMM-Newton, Swift, and NuS- 
TAR are still accumulating exposures on TeV blazars. For ex- 
ample, XMM-Newton has 9 Ms of exposure on TeV balazars 
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FIG. A.— The candidate rapid X-ray flare of LES 1101-232. The ObsID 
is 95387-02-07-00, and only PCU2 was in operation during the observation. 
The grey points and blue squares are the data of 16 s bins and 96 s bins, 
respectively. The red dashed curve is the least square model fitting to the 
96 s-bin light curve. 
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as of July 2017. A natural follow-up work is to search for ad- 
ditional rapid X-ray flares from TeV blazars in these archival 
databases. With the increased sample, we may study more 
details of the spectral variability of the flares, thus shedding 
light on their origin and particle acceleration processes. We 
can also study the correlations between flares and the proper- 
ties of the jets and central engine, as well as the reason for the 
small flaring timescales. 
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APPENDIX 
A. A CANDIDATE FAST FLARE OF 1ES 1101-232 


We found one fast X-ray flaring event of LES 1101-232, 
which is an HBL at a relatively high redshift (z = 0.186) com- 
pared with other HBLs. The full-band (2—20 keV) lightcurve 
in 16 s bins (Fig. shows elevated flux in the second half 
of the observation. When shown in 96 s bins, the lightcurve 
clearly manifests an almost complete flare profile, and can be 
fitted using the model of Eq.|1| The amplitude (F/F, = 1.93) 
of this flare is modest, and the rising timescale (7, = 60s) is 
almost as fast as PKS 2005-489. We only reported it as a 
candidate because the flux preceding the flare seems to be in 
continuous declining. The noisy data prevent us from per- 
forming data analysis in the same fashion as other sources 
(sub-band lightcurve fitting, joint spectral fitting, etc.). 


B. SPECTRAL VARIABILITY OF MRK 421 


The flare-only component of PKS 2005-489 shows an in- 
teresting spectral variability pattern, but the noisy data prevent 
us from making strong conclusions. The oblique “8” pattern 
in Fig. Dia is reminiscent of a flare of Mrk 421, which allows 
us to confirm our finding in another source using data from a 
different satellite. 
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FIG. B1.— An X-ray flare of Mrk 421 reported pp eee E 
Top panel: The dots and squares are soft-band and hard-band lightcurves 
respectively (cf. Fig. 6 of 2010), where the hard band has been 
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multiplied by a factor of 3. dashed lines are the constant components 
in two bands. In the bottom panel, the squares are the hardness ratios of the 
flare-only component. The smooth curves are model lightcurves (top panel) 


or model predicted hardness ratios (bottom panel). 
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FIG. B2.— Left: HR—flux diagram before the subtraction of the constant 
component (cf. Fig. 6 of (Garson tl [2010 The arrows indicate the di- 
rection of variation. Right: HR—flux diagram after the subtraction of constant 
component. Smooth curves are calculated from the model. The effective pho- 
ton index I’ on the y-axis of right hand side is obtained using the response files 
of XIS, Galactic absorbed power-law model and a range of assumed photon 
indices. 


We downloaded Suzaku data (ObsID 703043010) of 
Mrk 421 SESE E We only used data from the 
X-ray Imaging Spectrometer (XIS) onboard Suzaku. The data 
were reprocessed and screened using AEPIPELINE included 
in the Suzaku FTOOLS. In addition to the standard screen- 
ing criteria, we also required the cutoff rigidity to be larger 
than 6 GV/c following 010), We extracted 
lightcurves separately from the cleaned event files of XISO 
and XIS3 using XSELECT in initial 16 s bins. The source re- 
gion has an inner radius of 35 pixels and outer radius of 408 
pixels, while the background region is an annulus with an in- 
ner radius of 432 pixels and outer radius of 464 pixels. The 
lightcurves were extracted in two energy bands: 0.5-2 keV 
and 2—10 keV, which are referred to as the soft band and hard 
band, respectively. Lightcurves of the same energy band but 
different detectors (XISO and XIS3) were then merged. We 
rebinned the lightcurves in a size of 5760 s in accordance 


with the orbital period of Suzaku. Finally, the background 
lightcurves were subtracted to obtain the estimation of net 
count rates from the source. 

We clipped the lightcurves and kept only the segment of 
Flare 2 that is defined in{Garson et al.|(2010). We fitted the 
soft- and hard-band flares using a model that is analogous to 
Eq.|1} Only the constant component underlying the flare was 
replaced by a slowly-varying component using a linear func- 
tion Fk, = m+Slope xt. This component is in long-term de- 
cline in both bands. The fitted models were shown in the top 
panel of Fig. The error bars of the lightcurves are ex- 
ceedingly small, which renders the fitting statistically unac- 
ceptable and suggests ultimate incorrectness of the model due 
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to the existence of sub-structures. However, variation of the 
smooth model curves in the figure match reasonably well with 
observational data. 

The timescales are 7, = 9630 s and Ta = 9261 s in the soft 
band and 7, = 14756 s and Tq = 6796 s in the hard band, which 
is consistent with the pattern that hard band rises slowly and 
decays fast. The fitted long-term declining background was 
subtracted from the observed lightcurve of each band to leave 
out the flare-only component. We plotted hardness ratio vari- 
ation with time and flux in the bottom panel of Fig. and 
the right panel of Fig. Fig. shows the same pat- 
terns as panels (a)(b) of ba äi The effective photon indices 
(Ttf ~ 1.84) are also close to the photon indices of the flaring 
components that are listed in Table|3| 
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